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Abstract 

Heckler has recently argued that the Hawking radiation emitted from microscopic black holes 
has sufficiently strong interactions above a certain critical temperature that it forms a photosphere, 
analogous to that of the sun. In this case, the visible radiation is much cooler than the central 
temperature at the Schwarzschild radius, in contrast to the naive expectation for the observable 
spectrum. We investigate these ideas more quantitatively by solving the Boltzmann equation using 
the test particle method. Wc confirm that at least two kinds of photospheres may form: a quark- 
gluon plasma for black holes of mass Mbh £ 5 x 10 14 g and an electron-positron-photon plasma for 
Mbh ~ 2 x 10 12 g. The QCD photosphere extends from the black hole horizon to a distance of 0.2-4.0 
fm for 10 9 g < Mbh £ 5 x 10 14 g, at which point quarks and gluons with average energy of order 
Aqcd hadronize. The QED photosphere starts at a distance of approximately 700 black hole radii and 
dissipates at about 400 fm, where the average energy of the emitted electrons, positrons and photons 
is inversely proportional to the black hole temperature, and significantly higher than was found by 
Heckler. The consequences of these photospheres for the cosmic diffuse gamma ray and antiproton 
backgrounds are discussed: bounds on the black hole contribution to the density of the universe are 
slightly weakened. 
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1 Introduction 

It has long been known that black holes are 
not perfectly black, but emit nearly black-body 
radiation at a temperature Tbh = (87rGA/) _1 
due to quantum mechanical effects jlj (where G is 
Newton's gravitational constant). Although this 
Hawking radiation is negligible for astrophysically 
large black holes, it becomes sufficiently hot to 
be visible for masses M <, 10 15 GeV, correspond- 
ing to a BH that would be disappearing today, 
assuming it was present at the big bang. The 
present density of such BH's in the universe is in 
fact limited by observations of the diffuse gamma 
ray background coming from their accumulated 
radiation (for a review see §). Such limits arc 
calculated assuming Hawking's expression for the 
spectrum, 

dN _ a s {E)E 2 1 
dEdt 2ir exp(E/T BH )±l 

where the sign depends on whether the emitted 
particle is afermion (+) or a boson (— ), and <j s {E) 
is the absorption cross section for the emitted par- 
ticle, which depends on its spin s Moreover, 
previous estimates of the possibility of observing 
individual black holes, which explode in a burst 
of radiation as their masses approach to Planck 
mass, are also based on (|j), combined with calcu- 
lations of the mass spectrum of primordial black 
holes (PBH's) that could form during inflation Q 
or the QCD phase transition ||, and guesses as 
to how they might cluster. It is possible that 
extremely high energy neutrinos from exploding 
PBH's will be observable in the new generation of 
neutrino telescopes 0. 

Recently Heckler [g revived the possibility 
that the spectrum (|I|) need not hold far away from 
the Schwarschild radius, r# — 2GM, because the 
radiation might interact with itself in some re- 
gion, similar to photons diffusing inside the pho- 
tosphere of the sun. The idea had been previously 
dismissed jio], |), but in the framework of QED, 
Heckler identified bremsstrahlung and pair pro- 
duction as processes which could change this neg- 
ative conclusion. Both interactions cause a small 
initial number of high energy particles to fragment 
into many lower energy particles, giving a signifi- 
cant decrease in the average particle energy. Con- 
sequently, according to ref . || , a BH with T = 45 



GeV at th would appear to an observer more like 
a BH with T = m e (the electron mass), in terms 
of average energy, but much brighter in terms of 
absolute luminosity. 

As recognized by Heckler, it is not easy to 
give a highly quantitative treatment of the prob- 
lem because, unlike in the sun, the system is never 
so strongly coupled as to admit the approximation 
of local thermal equilibrium. The problem is that, 
whereas the sun has a huge chemical potential in 
gravitationally bound electrons and protons which 
are providing the large density, in the BH there 
are equal numbers of particles and antiparticles, 
which are at first freely streaming away from the 
horizon. The density at the horizon is ~ T 3 but 
quickly falls like 1/r 2 in the absence of particle 
production. In fact the density changes by a large 
factor within the mean free path of the particles. 
A fluid description of the plasma, although per- 
haps useful for roughly estimating its behavior, is 
not self-consistent. 

It thus seems worthwhile to investigate the 
evolution of the Hawking radiation plasma more 
quantitatively. The proper framework for doing 
so is the Boltzmann equation for the particle dis- 
tributions, 

(! + f -v)/(p,x,i)=C[/(p,x,i)] (2) 

using (Q) as a boundary condition at r = ru- 
in the following section we will discuss a gen- 
eral method for solving the Boltzmann equation 
which has been successfully used in the field of 
heavy ion collisions, and the adaptations of this 
method which we have made for the BH prob- 
lem. Crucial for this investigation is the collision 
term C[f]. Section 3 focuses on the scattering 
cross sections which go into C. In section 4 we 
present our results for the detailed properties of 
the QED and QCD photospheres. The ramifi- 
cations for the most relevant observable particle 
backgrounds, namely gamma rays and antipro- 
tons, are worked out in section 5. We summarize 
our results in the final section. 
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2 Test particle method for 
solving Boltzmann equa- 
tion 

The full Boltzmann equation is an integro- 
diffcrential equation which is difficult to solve ex- 
actly. For black holes we can make several impor- 
tant simplifications. (1) We consider nonrotating 
black holes so that the distribution functions have 
spherical symmetry. (2) We confine our attention 
to BH's whose lifetime is still long compared to 
the diffusion time of particles in the putative pho- 
tosphere; thus the distribution functions are ap- 
proximately time-independent. For example, the 
lifetime of a black hole with M — 10 11 g and 
T = 1000 GcV is of order 10 5 s, whereas the par- 
ticle diffusion time in its photosphere is ~ 10~ 21 
s. We estimate that only for black holes of mass 
<, 10 7 g will the lifetime become comparable to 
the diffusion time. 

Under the above assumptions, the Boltzmann 
equation takes the simpler form 

^-f(p, Pt ,r) = -C[f}, (3) 
or v r 

where p = |p|, p t is the magnitude of the compo- 
nent of momentum transverse to the radial direc- 
tion, and v r is the radial velocity. Because of the 
spherical symmetry there is no dependence on the 
azimuthal momentum component. 

Even with these simplifications, eq. (|^) is 
still prohibitively difficult to solve in the most 
naive numerical way, namely discretizing momen- 
tum space and evolving the distribution defined 
on this momentum lattice forward in radius. At 
each point in momentum space a multidimen- 
sional phase space integral is required to evalu- 
ate the collision term, which is computationally 
costly. Fortunately, this problem has been already 
surmounted in other situations, namely heavy ion 
collisions [|ll[ There one wants to track the 
distribution functions of nucleons in the two nu- 
clei as they pass through each other and undergo 
collisions. 

The idea is essentially to follow the classi- 
cal evolution of each particle, allowing for possi- 
ble scatterings by using the differential scattering 
cross section as a probability distribution. How- 
ever, the number of nucleons in even a heavy nu- 



cleus is so small that one must do this many times 
to obtain distribution functions that are not dom- 
inated by statistical fluctuations. Alternatively, 
one can obtain the required statistics by repre- 
senting each real particle by an arbitrarily large 
number of test particles, which is equivalent to 
but simpler than simulating the heavy ion colli- 
sion many times. One must only be sure to avoid 
unphysical collisions by prohibiting any two test 
particles that represent the same actual particle 
from scattering off each other. However, in prac- 
tice, it is easier to allow scatterings to occur be- 
tween all test particles, while simultaneously re- 
ducing the cross section by a factor N, equal to 
the ratio between the number of test particles and 
the number of actual particles. 

A difference between heavy ion collisions and 
black hole radiation is that in the former, one 
is evolving the distributions in time, whereas we 
want to evolve them in space. However, our ver- 
sion of the Boltzmann equation is mathematically 
equivalent to one with time evolution but spatial 
homogeneity, 

|/(p, *)=£'[/], (4) 

where /(p, t) is known at some initial time. Since 
this equation is evidently amenable to solution by 
the test particle method, and eq. ([}]) is equivalent 
to it just by renaming t — > r and C'[f] — > C[f]/v r , 
we can also apply the same method to solve (||). 

2.1 Free evolution 

Let us first consider how the method works in the 
absence of collisions, so the particles are propa- 
gating freely. At the initial surface of the horizon, 
it is assumed that the radiation is distributed ac- 
cording to eq. (Q) in all directions not pointing 
back below the horizon, that is, for polar angles 
with cos# > 0. However, at a larger radius r, if 
there have been no collisions then particles can 
only have come from a cone pointing back to and 
subtended by the horizon. The directions of the 
possible momenta of the particles are restricted to 
lie in the outward half of this cone, whose opening 
angle is given by sm9 — ru/r. As r increases, the 
momenta become increasingly concentrated in the 
radial direction. 
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Numerically the evolution in the noninteract- 
ing case is thus trivial. A particle whose momen- 
tum is p at rji will have the same momentum at 
r' > r. However, the coordinates of the momenta 
will change: for example, if p is purely transverse 
at r, it will develop a radial component at r' . Let 
9 be the polar angle from the radial direction at r 
and 8' that at r' (see Fig. 0). Then using the law 
of sines, 
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2.2 Including interactions 

To account for the interactions of the particles 
with each other one must first choose an appro- 
priate step size Ar for evolving the distribution. 
The natural choice is the mean free path, A, or 
some fraction thereof. If the interaction cross sec- 
tion is (7, then we can define 



A x (r) = n(r)(av re i), 



(9) 



where n is the density and v re i is the relative ve- 
locity of the two interacting particles (to be dis- 
cussed in more detail below, in section 3). How- 
ever, the time-independent Boltzmann equation 
we are solving is mathematically equivalent to the 
usual time-dependent version only after dividing 
by the radial velocity. Thus our collision term is 
related to the usual one by a factor of 1 jv r , which 
modifies the definition of the mean free path to 



Figure 1: Momentum in free radial evolution 



\-\r)=n{r)(n^ 

V r 



(10) 



For very large r' the momentum becomes es- 
sentially radial, 8' = 0, and a momentum that is 
initially radial remains so. The distribution func- 
tion will be given approximately by 



/(p'X) 



f(p,r H ), 0<sin8' <r H /j 
0, sin0' > r H /r 



(6) 



where the new radial and perpendicular compo- 
nents arc related to the old ones by 



Pi 



—P± and p' r = \J p 2 - p' ± 2 . (7) 



If / is isotropic at the horizon, then for large dis- 
tances one finds that this gives a density decreas- 
ing like 1/r , as required by conservation of par- 
ticle flux: 



n{r) 



70/(P-O = (^) 2 nM. (8) 



Since the particles are moving on straight lines, 
the step size for evolving the distributions is irrel- 
evant in the noninteracting case. 



We take v r in this formula to be the minimum of 
the radial velocities of the two interacting parti- 
cles. In practice this factor is not very important 
because the bulk flow of the particles away from 
the black hole tends to cause v r to approach unity. 

The average is taken using the distribution 
functions for the incoming particles, 



Vrel 



d 3 p 1 d 3 p 2 
(2tt) 6 



Vrel 



f(Px,r)f(p2,r)a 



(11) 

In the test particle method, there is no need to 
do any integrals per se. Since the ensemble of 
test particles is already distributed according to 
/, (o-v r ei/vr} is simply given by an unweighted av- 
erage of av re i/v r over the ensemble. We randomly 
choose pairs of particles to perform this average. 

If A(r) was constant, one could choose the 
step size for the evolution to be simply Ar = A. 
Since A is the average distance particles go be- 
tween collisions, the interactions could be simu- 
lated by demanding that over the distance Ar = 
A, each particle participates in a single collision. 
However, in the black hole problem A(r) can in- 
crease significantly on a distance scale of A be- 
cause the density is decreasing like 1/r 2 . To deal 
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with this, one must choose the step size to be 
smaller, 



Ar < 



A 



dX/dr' 



(12) 



, electrons 
photons. 



Then the density at the horizon is 



(15) 



Then A, and hence the interaction rate, is guaran- 
teed to be approximately constant over Ar. Over 
this distance, only a fraction F of all the particles 
will undergo scattering, 



(13) 



The general procedure for including interac- 
tions is therefore clear: (1) At each r, compute A 
using the known distributions /(p, r). (2) Given 
A, choose a step size Ar in accordance with (|l2|). 
(3) Choose a subset of particles from the ensem- 
ble such that a fraction F = Ar/A of the total 
ensemble will undergo collisions and change their 
energies and momenta accordingly. To account for 
the l/v r factor in the effective cross section, we ar- 
range the ensemble in order of increasing v r and 
let that fraction of particles F with the smallest v r 
participate in the interactions. After these steps 
are carried out, the distribution function is known 
at r + Ar and the procedure can be repeated. 

The particle density n(r) which goes into ( [To| ) 
is calculated analytically and therefore does not 
depend on the number of the test particles. To 
find the density n(r) we first compute what it 
would be in the absence of particle production: 



n (r) = n/j-|, 



(14) 



where the subscript means that this density 
is before the particle creation process is taken 
into account. The radius of the horizon is 
r/j = 1/(4:1tTbh)- To find the density at the hori- 
zon, we note that if the BH absorption cross sec- 
tion a(E) was a constant, then the BH would be 
a perfect black body, and the density of radia- 
tion would be thermal. However a(E) vanishes 
as E — > and only reaches its geometric optics 
value of do = 27nG 2 M 2 in the limit as E — + oo. 
Therefore rih is reduced from its thermal value by 
a factor of T s = J d 3 pa(E)f±(p)/ J d 3 pa f±(p). 
This has been computed in ref. M to be 



2§jT/ C(3) Tg H , electrons or positrons 

"•h = { 

^nC(3)T 3 H , photons, 

(16) 

where £(3) = 1.20206 (Rieman zeta function). 

To account for particle creation we use the 
test particles to find the fraction of new fermions 
and bosons created at each step. Let Nf^r) be 
the number of electrons [photons] in the shell of 
width Ar at radius r. We will define Pyrw(r) as: 



Then, using (|lj) through ( [Tt| ) we obtain 



(17) 



n(r) 



C(3) T BH 
7r 2 (47r) 2 r 2 




(18) 



In this way, particle densities can be computed at 
any step by keeping track of the relative increase 
in particle number, P/rw- Later, we shall also re- 
fer to the ratio of all particles at r versus at the 
horizon, 



P(r) 



n b (r) + rt/(r) 
n b {r h ) + n f (r h ) 
P b {r) +4.1 7P/(r) 
1+4.17 



(19) 



The factor 4.17 comes from computing STf/ZT},. 

To generalize the previous results to quarks 
and gluons is straightforward: one must multiply 
the photon density by a factor of 8 to convert to 
gluons, and the electron density by a factor 3n/ 
to get nt flavors of quarks. 

A final issue concerns the number of test par- 
ticles used to represent the ensemble, versus the 
actual number of particles coming from the BH. 
For example, within the first radial increment Ar 
near the horizon, the actual number of particles 
AN — Aitr 2 H Ar n{rH) might be too small a num- 
ber to generate a smooth distribution function. 
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We would prefer to represent these particles with 
some much larger number N t of test particles. In 
the application to heavy ion collisions it is nec- 
essary to reduce the cross section by the factor 
AN/N t to avoid overestimating the number of col- 
lisions. This is because in the latter situation, in- 
dividual nucleon positions are kept track of, and 
two particles are only allowed to collide if they 
come within a distance b c = \fo~pK of each other. 
If the number of particles is artificially increased 
while proportionally decreasing the cross section, 
the physically meaningful mean free path will re- 
main constant. In our case, since we do not follow 
the spatial trajectories of particles, but instead 
compute the physical mean free path, it is con- 
sistent to allow all the test particles to interact 
over a distance A, regardless of how large the ra- 
tio N t /AN is. Thus there is no need to reduce a 
proportionally to the number of test particles in 
the present problem. 

3 Scattering cross section 

The most important processes contributing to 
photosphere formation around the BH are 
bremsstrahlung (ee — > ee~/) and photon-electron 
pair production (ey — > ee~e + ), whose dominant 
Feynman diagrams are shown in Fig. 0. Both 




Figure 2: Feynman diagrams for the dominant 
contributions to (a) bremsstrahlung ; (b) pair pro- 
duction. 

have cross sections of 0(a 3 ), which at first sight 
might make them seem less important than the 
leading 0(a 2 ) processes like Compton scattering. 
However the latter are elastic and change the 
distributions only by randomizing the momenta, 
whereas the former are inelastic and increase the 
number of particles while reducing their average 
energies. Thus the elastic scattering processes ap- 
pear to be less relevant for determining the prop- 
erties of the photosphere even though they are 
faster than the inelastic ones. Before going on to 
a detailed description of brchmsstrahlung, let us 
examine this issue more carefully. 



3.1 Elastic scattering processes 

One could argue that elastic interactions like 
Compton scattering (cy — > ej) or electron- 
electron scattering might change the size of the 
photosphere somewhat, but they will not affect 
the most important observable property, which is 
the energy spectrum of emitted particles. It is 
the 2^3 body interactions like bremsstrahlung 
which principally distinguish the photosphere 
from freely streaming particles. However there is 
one sense in which elastic scatterings may be im- 
portant: by randomizing the particle momenta, 
they might postpone the tendency for the relative 
velocities of particles to approach zero, due to the 
bulk, outward radial motion. This in turn could 
keep the bremsstrahlung mean free path, eq. ([lO]) , 
small out to larger radii, enlarging the size of the 
photosphere. 

Suppose the mean free path for elastic pro- 
cesses is given by A e (r) at a distance r from the 
black hole. Their momentum-randomizing effects 
will be important for the photosphere only in re- 
gions where A e < r. This is because the length 
scale over which random momenta become in- 
creasingly radial at a given distance r is r itself. 
If Ae(r-) exceeds r, then scatterings start to lose in 
the competition against geometry. 

The relativistic limit of the Compton cross 
section, in the center of mass frame, is 

2?ra 2 , E 
a c = — In . 20 

m e E m e 

We can then estimate A e as 

A" 1 ^) S a c n{r) = a c n(r h )P(r), (21) 

with the particle density given by (|l^) and the 
particle production factor by (|l9|). Below, we 
will show that when elastic scattering is ignored, 
the QED photosphere ends at a distance of rt ~ 
m" 1 , and for large Tbh, P{ r f) ~ bT BH w ^ n 
b = 4.5 x 1CP 5 GeV~ 2 . Furthermore, the aver- 
age particle energy at the horizon will be given by 
its value at the horizon, ~ 3Tb h, divided by P(r). 
To see whether our neglect of elastic scattering is 
consistent, we should compute rf/X c = (m e \ c )~ l , 
and ask whether it ever exceeds O(l). Using the 
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above results we get 




Therefore the effects of elastic scattering should 
only become important when the BH tempera- 
ture starts to exceed ~ 5 TeV. We will discuss 
our numerical investigation of this regime below, 
although most of our work focuses on BH temper- 
atures below 1 TeV. 

3.2 Bremsstrahlung cross section 

The relativistic differential cross section for 
bremsstrahlung in the center of mass frame is 

du Eu \i K ' E ) 

x In K - >- - - , 23 

V mjuj 2 / 

where K = c = 1, r e = a/m e , E is the initial 
energy of each electron, and u> is the energy of 
the emitted photon. It diverges for a; — » 0, but 
higher order corrections essentially impose an in- 
frared cutoff ||. 

The form of Eq. (|23|) implies that the proba- 
bility to emit a photon diverges as its energy goes 
to zero. On the other hand, emission of a zero- 
energy photon has no effect on the electron which 
emits it. A convenient way of rendering the cross 
section finite, while at the same time accounting 
for the photons which carry away significant en- 
ergy, is to use the energy-averaged total cross sec- 

tion §, m 

^Je^)^* 8 ^^- (24) 

The cross section for photon-electron pair produc- 
tion shows the same functional dependence in the 
extreme relativistic limit and we therefore use the 
same estimate ( pij) for both interactions. 

An improvement on the present treat- 
ment would be to use the actual differential 
bremsstrahlung cross section as a distribution 



which would produce low energy photons with 
higher probability. One would have to impose 
an infrared cutoff on the emitted photon energy 
and show that no meaningful physical properties 
of the photosphere depend on this cutoff. Here 
we have taken the simpler approach of approxi- 
mating da/dCl — a /Ait, and choosing the energies 
and directions of the final state particles at ran- 
dom in the center of mass system, subject to the 
constraint of conservation of four-momentum. 

3.3 Thermal Mass 

Bremsstrahlung proceses do not occur in vacuum 
but in a background plasma of almost radially 
propagating particles. These background parti- 
cles suppress the bremsstrahlung cross section if 
they are sufficiently dense (the LPM effect p(i[). 
A simplified means of accounting for this is to re- 
place the vacuum electron mass, m e , by its ther- 
mal value, 

m th = y/ml + ml(T). (25) 

This procedure, although not exact, gives the cor- 
rect position of the pole of the electron propaga- 
tor in the energy-momentum plane to an accu- 
racy of 10% pi. In a gauge theory with coupling 
constant g, in a thermal background, the plasma 
mass is given by m p = gT^C(R)/8, where C(R) 
is the quadratic Casimir invariant for the fermion 
transforming in the representation R of the gauge 
group. However, even if the background is not 
in thermal equilibrium, we find by computing the 
electron self-energy in the plasma that the formula 
for m p which was originally derived in the thermal 
case is still correct: 




(ff(p) + fb(p)), (26) 



where //(6)(p) is the distribution function for a 
single polarization of the background fermions 
(gauge bosons), (i.e., f f ( b) {p) = {e E/T ± I) -1 in 
the thermal case). For a collection of test parti- 
cles which represents these distribution functions, 
the integral may be approximated by taking the 
average of 1/p over all electrons and photons, re- 
spectively (quarks and gluons in the case of QCD). 
Counting polarizations, we obtain 



G 



Ana (n e (p e x ) +n 7 (p 7 1 )) , 



QED 
QCD. 



(27) 



Near the BH horizon where the particle den- 
sities are high, the plasma mass dominates over 
the vacuum mass, and leads to a significant re- 
duction in the bremsstrahlung cross section, since 
a oc m^ 2 . Failure to take this into account would 
give a misleading estimate of the onset of the pho- 
tosphere in the case of QED. 

Several aspects of the implementation of ther- 
mal mass in the test particle model should be men- 
tioned. First, in eq. ( p7| ) n e and n 7 are the actual 
electron and photon densities calculated analyt- 
ically, whereas p e and p 7 arc kinetic energies of 
the test particles, assumed to be relativistic. Sec- 
ondly, we recalculate the fermion thermal mass 
[using eq. (^||)] at each step in radius. Because 
m t h changes, if the fermion momentum was held 
fixed, its total energy E — y m th + Pe ; would 
change. Hence, we assume that the thermal mass 
correction is done at the expense of momentum in 
such a way as to keep its total energy E constant. 

(The reader may wonder whether energy 
should indeed be conserved, due to gravitational 
redshift. In fact all significant redshifting of the 
outgoing particle energies occurs within the first 
few Schwarzschild radii of the horizon, where the 
classical concept of particles is not yet well de- 
fined, due in part to the deBroglie wavelength be- 
ing of the same order as the Schwarzschild radius. 
Moreover, the Hawking spectrum (Jl]) refers to en- 
ergies as observed far from the BH. If we were 
going to try to include redshift effects, we should 
blueshift the initial particle energies accordingly, 
so that asymptotically they have the usual dis- 
tribution. This would only affect our simulation 
very close to the horizon, which does not appear 
to be an important region as regards the observ- 
able features of the photosphere. Hence it is sim- 
pler and seemingly a good approximation to ig- 
nore redshift.) 

As in ref. ||, we can estimate the ther- 
mal mass of an electron in the QED photo- 
sphere using n(r) = T BH (3/2)^ /(7r 2 (47r)V) 
for the densities of both photons and electrons, 
and 3T Bff /(3/2)^( r ) for average particle energy 
at radius r. Here, J\f(r) is the number of 
bremsstrahlung events an average particle under- 
goes between the horizon and radius r. Since in 



each such scattering one new particle is produced, 
(3/2) M ^ is the factor by which the total num- 
ber of particles has increased. In the test-particle 
model, since the total number of particles N(r) 
is known at each step, this increase is calculated 
directly as N(r)/N(r h ). Then, 



m t h 



a 1 / N(r) 



ir(4ir) 2 r 2 \N(r h ) 



(28) 
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Figure 3: Thermal mass (and the estimate (28)) 
in typical (a) QED and (b) QCD photospheres. 

Fig. |^ shows how the thermal mass depends 
on radial distance from the BH for two illustra- 
tive cases: the QED photosphere of a BH with 
T = 1000 GeV, and the QCD photosphere for a 
T = 10 GeV black hole. We see that in the QED 
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plasma the factor of 1/r 2 is dominant over any 
growth in density as the radius increases, which 
causes the thermal mass to decrease monotoni- 
cally. On the other hand, in the QCD plasma, 
the particle production factor N(r) /iV(r^) is more 
influential (due to the running of the coupling 
constant with energy, in particular) and just be- 
fore hadronization starts, at average particle ener- 
gies (E) near the QCD scale Aqcd, the thermal 
mass may again increase. Our perturbative for- 
mulas should not be trusted quantitatively in this 



regime, however. Fig. 
lytic approximation (t 



3 a also shows that the ana- 
3) rath and the numerically 



computed values are in reasonably good agree- 
ment. 

3.4 Interaction Distance and Rela- 
tive Velocity 

The distance at which particles can interact via 
bremsstrahlung and pair production is an impor- 
tant parameter, since we have to decide which test 
particles may interact with one another and which 
may not. In the heavy ion collision problem, those 
test particles were allowed to interact which were 
within a critical distance b c = y 'a /tt from each 
other. Notice that in that case, the relative veloc- 
ity of the two particles was typically large because 
the interacting nucleons were not moving parallel 
to each other. But in the present case of BH ra- 
diation, two particles which are nearest neighbors 
typically are moving in almost the same direction, 
i.e., radially, leading to a suppressed relative ve- 
locity v re i. Such particles have a small probabil- 
ity of interacting, since they make a small con- 
tribution to the inverse mean free path, eq. (||). 
The dominant interactions involve pairs of parti- 
cles with large relative velocity. Near the horizon 
of the BH, such pairs would consist of particles on 
opposite sides of the BH, moving away from each 
other. 

Although classically it is somewhat counter- 
intuitive to have interactions between particles 
that have already passed each other, as it were, 
there is no reason to exclude them so long as the 
particles are still within the range of the inter- 
action. For bremsstrahlung it has been shown 
[fTcf that at relativistic energies, small momentum 
transfers of the order k ~ m e contribute the bulk 
of the total cross section. Therefore the distance 
at which particles can interact via bremsstrahlung 



(and pair production) is of order m" 1 in the vac- 
uum. In the plasma, accounting for the LPM ef- 
fect, the cutoff instead becomes the thermal mass 

Therefore in the radial interval < r < m th 
(assuming the BH is microscopic so that the pho- 
tosphere can indeed form), a given particle in the 
plasma is always able to interact with some other 
particle, with a large relative velocity, v re i ;$ 2. 
It is not necessary to assume particles interact 
only with their nearest neighbors, or to keep track 
of the exact spatial trajectories of the particles, 
which simplifies our task of evolving the ensemble 
of test particles. 

Of course the above procedure no longer 
works at radii larger than m7 h . At that point 
we adopt the procedure of randomly pairing par- 
ticles, assuming the pairs represent nearest neigh- 
bors separated by the average interparticle dis- 
tance n(r) 1//3 , and computing their actual rela- 
tive velocity, given by the formula 



Vrei = 



y(px ■V2) 2 - m\m\ 
E\E<i 



Pi • P2 

PlP2 




(29) 



whose second form is valid in the relativistic 
limit .^| In this large-radius regime we rely upon 
the randomization of velocities provided by the 
scatterings themselves to keep the momenta from 
becoming purely radial, which otherwise causes 
Vrei to tend toward zero. Nevertheless the in- 
evitable radialization of momenta as r increases 
quickly overcomes the randomization due to scat- 
terings, causing the particles to become freely 
propagating and marking the end of the photo- 
sphere. This behavior is shown in Fig. [|. Al- 
though our matching between the r < m~^ h x and 
r > m"^ regions is somewhat crude, we be- 
lieve it captures the essential physics, and that a 
more accurate treatment, using detailed informa- 
tion about each particle's trajectory, would only 
change the estimate of the photosphere's outer ra- 
dius by a factor of order unity. 

2 This definition of v re i arises from comparing the de- 
tailed form of the collision term in the Boltzmann equation 
with that of the scattering cross section. It is important 
to notice that the scattering angle between pi and p2 is 
evaluated in the rest frame of the BH. This is the same 
frame in which the Boltzmann equation is most naturally 
formulated for the present problem. 
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Figure 4: Maximum distance of interaction, m7 h , 
in a typical QED photosphere (Tb# = 1000 GeV), 
indicating the rapid dissipation of the photosphere 
at distances r > m^ h . 



Heckler The idea is to dehne a quantity Af(r) 
denoting the number of scatterings an average 
particle undergoes between the horizon and some 
radius r > T%. The criterion for having a photo- 
sphere is taken to be that on average every particle 
undergoes a collision at least once between leav- 
ing the horizon and escaping to infinity, in other 
words that 



lim Mir) > 1. 



(30) 



T c is then the temperature of a black hole for 
which this limit is exactly unity. The critical tem- 
perature we obtain is T c ~ 50 GeV, whereas the 
result following from the approximate method of 
ref. H is 45.2 GeV. Our determination of Af(r) 
is shown for several black hole temperatures in 
Fig. | 
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4 Numerical Results 
QED photosphere 

This section will summarize the results of test par- 
ticle simulation of black hole emission when only 
electrons, positrons and photons are taken into 
consideration. This restriction is appropriate for 
black holes with Tbh < ^-qcd- We will also treat 
higher temperatures in this context for the sake 
of understanding, deferring study of the effects of 
quarks and gluons until the next subsection. 
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4.1 Photosphere Formation 

We confirm the formation of a photosphere for 
black hole temperatures Tbh ~ 100 GeV and 
higher, as was originally suggested in ||. The 
photosphere first appears at a radius r ~ 10 4 r^ 
(where r^ — 1/AttTbh is the Schwarzschild ra- 
dius), and is characterized by a region of in- 
tense collisions terminating at a distance of rf ~ 
10 7 — 10 8 rh from the black hole. The effective- 
ness of the collisions is demonstrated by a very 
slow increase, or even decrease, of the mean free 
path in the photosphere compared to the interior 
and exterior regions. 

The photosphere forms only for black holes 
above a certain critical temperature T c . We use 
the same definition of T c as was introduced by 



Figure 5: The number of scatterings an average 
particle undergoes between the horizon and radius 
r: N{r). 



4.2 Photosphere Parameters 

Here we will present the photosphere parameters 
obtained from simulations for approximately 30 
different black hole temperatures. These parame- 
ters include the radii of the inner and outer sur- 
faces, Ti and rf respectively, the total particle pro- 
duction factor P and the average energy of the 
particles emerging from the photosphere E / . The 
latter is relevant because it is the average energy 
of particles that may eventually reach a distant 
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observer. We will point out significant discrep- 
ancies between these results and the fluid model 
used in Q , and derive empirical formulas from our 
simulations showing the dependence of fj, Tf and 
Ef on the BH temperature, Tbh- 

Inner radius 

The radius of the inner surface of the photosphere 
(r.;) is defined by N{r,i) = 1, i.e., the radius by 
which on average every particle originating at the 



dieted (fj ~ lO 9 ^) in M using a fluid model for 
the interacting particles. 

Edge Radius 

The outer radius of the photosphere, r/, can be 
defined as the distance after which particles effec- 
tively free stream away without significant inter- 
actions with one another. The mean free path A 
quickly begins to diverge at this point. If A^ is 
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Figure 6: Radii of inner QED photosphere surface 
for different black hole temperatures (black dots). 
Solid line represents inverse regression over these 
values, see Eq. (f?l|). 

horizon has undergone one collision. The values 
of r*j in units of 1 GeV" 1 ~ 0.197 fm are plotted 
in Fig. H as a function of black hole temperature. 
As one can see from the graph, T{ decreases with 
the temperature and is closely fitted by 

n = -=—, k = (6.446±0.003) x 10~ 4 . (31) 

KJ-BH 

We know that the radius of the Schwarzschild 
horizon is also inversely proportional to the black 
hole temperature Th = 1/AttTbh, so that 

n = —r h ~2xl0 4 r h . (32) 

K 

By this criterion, the photosphere starts to de- 
velop much closer to the black hole than was pre- 



Figurc 7: Radii of outer QED photosphere surface 
for different black hole temperatures. The solid 
line shows r — l/m e , the radius beyond which 
particles can no longer interact with each other 
by bremsstrahlung . 

the mean free path at step i then the condition 
Ai+i/Aj ^> Xi/Xi-i is a convenient criterion to 
define the end the photosphere, and rj is just the 
radius the particles have reached by step i. This 
value is easily found in practice because in the 
next step after rj, A is usually several orders of 
magnitude larger than its previous value, imply- 
ing that the particles have become virtually free. 
Because A(r) is diverging so quickly near rj, the 
latter is inherently difficult to determine precisely, 
even when the steps size is made smaller near the 
edge of the photosphere. Great precision in r/ is 
not essential, however, because it is not an observ- 
able quantity. 

The values of r/ remain in the same range 
(1700-2200 GeV- 1 ) for all black hole tempera- 
tures shown and are consistent with being inde- 
pendent of Tbh at the horizon. Only the statis- 



10 



tical dispersion around the mean value decreases 
with rising temperature (Fig. fj]) . The mean value 
is approximately rf ~ 2000 GeV -1 , which is close 
to l/m e . This is easy to understand in light of our 
previous discussion (section 3.3 and fig. 3) of the 
range of the bremsstrahlung interaction, whose 
maximum value is of order l/m e . Because the 
trajectories of particles within the range of inter- 
action rapidly become parallel at radii beyond this 
cutoff, the photosphere quickly dissipates. 

Particle Production 

Another useful parameter for characterizing the 
photosphere is the total particle production fac- 
tor, given by P — N(rf)/N(rh), where N(r) is 
the number of test particles at radius r. P can 
be used to quantify the probability of interac- 
tions inside the photosphere, since at each col- 
lision N — ► N + 1. Figure || shows P as a function 
of black hole temperature. We find that it can be 
represented by a quadratic fit (solid line): 



where 



b = 



P{T B h) ^ clTbh + bT BH , 

0.026 GeV" 1 , T BH < 2 TeV; 

0.226 GeV" 1 , T BH > 2 TeV; 
4.5 x 10~ 5 GeV" 2 , T BH < 2 TeV; 

4.2 x 10" 5 GeV" 2 , T BH > 2 TeV. 



(33) 



(34) 



The dominant quadratic term implies that the 
particle number density in the photosphere in- 
creases rapidly in the last stages of black hole 
evaporation, as T approaches the Planck scale. 
The slightly different dependences for T BB > 2 
TeV and T BB < 2 TeV can be understood by 
the argument of section 3.1: this is the tempera- 
ture where randomizing effects of elastic scatter- 
ing become relevant. However, the effect is not 
dramatic. 

Average Final Energy 

From an observational point of view, the reduction 
in average energy of emitted particles is one of the 
most relevant consequences of the photosphere. 
At the horizon, Ei is approximately 3.1Tbh- But 





140 






o 


120 


o 




3 




sz 


100 


o 




uct 


80 


"D 




o 

1— 




Q. 


60 


CD 




O 






40 


CO 




Q_ 






20 








Results 

Quadratic regression 



200 400 600 800 1000 1200 1400 1600 

Black hole temperature, GeV 



4000 



£ 3500 
O 

™ 3000 
c 

■° 2500 
O 

"g 2000 
i_ 
Q. 

0) 1500 
o 

« 1000 
co 
CL 

500 



1000 2000 3000 4000 5000 6000 7000 8000 

Black hole temperature, GeV 



Figure 8: Particle production factor N(rf)/N(rh) 
versus T BH for (a) Tbh < 1.5 TeV and (b) T BH < 
8 TeV. Solid line is the quadratic fit to the results 
[eq. ©]. 



the photosphere can reduce this number dramati- 
cally, so that a distant observer sees a much softer 
spectrum. Our results for the Tg^-dependence of 
the average energy at the edge of the photosphere, 
Ef, are displayed in Figure ||. These values are 
significantly higher than those (~ 1 MeV) found 
in ref . || , and the temperature-dependence is also 
very different. This apparently stems from the use 
of a fluid description in j^] which is not really ap- 
plicable. 

We can compute E / from Ei and the particle 
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T t 




60 GeV 
300 GeV 
1000 GeV 


1.5 x 10 3 fm 
30 fm 
7.9 fm 


8.6 x 10 3 fm 
9.6 x 10 3 fm 
9.6 x 10 3 fm 


97.0 GeV 
79.5 GeV 
44.2 GeV 



Table 1: QED photosphere parameters (temperature at horizon, inner and outer photosphere radii, 
and average particle energy at the outer edge of the photosphere) for several black hole temperatures 
obtained in test particle simulations. 
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production factor, Eq. (|33| 



Ei = (2.7^ + 3.15^)T BH 
P 

3.1 



aT, 



BH 



bT 2 BH 



b + cT, 



BH 



(35) 



where rib, n f and n are the respective densities 
of bosons, fermions, and all particles at the hori- 
zon. In the limit of high black hole temperatures 
(a < bT BH «*■ T BH > 10 4 GeV): 
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(36) 
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Figure 9: Ef, average particle energy at the outer 
edge of the QED photosphere, for (a) Tbh < 1.5 
TeV and (b) T BH < 8 TeV. 



This is a remarkable result since it predicts that 
for black holes of temperature Tbh ^100 GeV 
(M < 10 12 g) the effective temperature of emit- 
ted particles (that which would be inferred by ob- 
servers far away from the hole) is actually lower, 
the higher is their temperature at the horizon. For 
an individual black hole, which is losing mass and 
hence becomes hotter, this means that its appar- 
ent temperature goes down. The word "temper- 
ature" is used loosely here: the spectrum is non- 
thermal, with a higher luminosity at low frequen- 
cies than that of a blackbody (see fig. 10 below). 

Eq. (|36| ) was derived on the assumption 
that the black hole horizon temperature does not 
change much on the time scale of particle diffusion 
from the horizon to the outer edge of the photo- 
sphere. Nevertheless, we expect the qualitative 
picture to be the same even for black holes tem- 
peratures above this limit of validity. It should 
also be kept in mind that we are discussing only 
electromagnetically interacting particles so far. 
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The behavior of the QCD photosphere (below) is 
quite different. 

The time development of an individual BH, 
based on the above results, can now be described. 
As far as only QED emission is concerned, a pho- 
tosphere starts to develop around the black hole 
when it evaporates to a mass of M ~ 5 x 10 12 g. 
At this point the average energy of emitted par- 
ticles, instead of going up inverse proportionally 
to the BH mass, levels off and begins to decrease. 
On the other hand, the total luminosity increases, 
and the spectrum becomes softer than that of a 
black body. The outer edge of the photosphere 
remains at a radius of 400 fm ~ m~ . How- 
ever, its inner radius, ~ 700 r^, shrinks with 
the Schwarzschild horizon rv Eventually, if the 
steady-state assumption remains valid at these 
temperatures, the edge will cool to Ef ~ m e , 
when the positrons and electrons annihilate and 
no further cooling occurs. At this point, however, 
the black hole is within 10~ 10 s of its total evap- 
oration. 

To give some idea of how the results of the 
test particle method differ from the estimates in 
where a nonperfect fluid model was used, we 
tabulate the relevant quantities in table |l|, for 
three different BH temperatures. 
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4.3 Inside the Photosphere 

A more detailed picture of the plasma can be seen 
from the particle momenta distributions at differ- 
ent radii inside the photosphere. Since no signifi- 
cant interactions start before r ~ 10 4 rh, particles 
move increasingly radially until the inner bound- 
ary of the photosphere. By this point they have on 
average undergone one interaction and their mo- 
menta begin to get randomized. Close to the end 
of the photosphere, the mean free path increases 
the radial components of the momenta again start 
to increase. These developments can be seen in 
the distributions of the transverse momenta. In 
fig. [lO] we show the distributions of the fraction 
of momentum which is transverse to the radial 
direction, pt/p — Pt/y/Pt + Pr> f° r several radii 
inside a sample photosphere. These figures show 
the overall growth in particle density as well as 
the shift to larger p t as one goes from the inner 
boundary of the photosphere toward its interior. 




0.0 0.1 0.2 0.3 



0.4 0.5 0.6 0.7 O.i 
P/P 



Close to the outer edge 




0.0 0.1 0.2 0.3 



0.4 0.5 0.6 0.7 
P/P 



Figure 10: Distribution of pt/p for several radial 
shells in the QED photosphere of a BH with T = 
1000 GeV. 
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4.4 Spectra 



Finally, we examine the particle energy spectra 
at the black hole horizon and at the photosphere 
edge (Fig. |ll|). The shift toward lower energies is 
the most significant difference between the origi- 
nal distribution at the horizon and those within 
the photosphere. The two are shown together in 
Figure [j^, where it can be seen that the soften- 
ing in energy is accompanied by an increase in 
the number of particles, as is required by energy 
conservation. 

The QED photosphere by itself serves as a 
kind of toy model for realistic black holes, which 
are also emitting quarks and gluons at the tem- 
peratures we are considering. We now turn our 
attention to the QCD photosphere. 
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Figure 11: Energy spectrum of a 1000 GeV 
(10 11 g) black hole at the Schwarzschild horizon 
and near the edge of the QED photosphere. 
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Figure 12: Same as figure 11, but with the two 
spectra superposed. 

QCD photosphere 

A black hole whose temperature exceeds Aqcd ~ 
200 MeV, thus having a mass M BH < 5 x 10 14 
g, emits quarks and gluons which, as proposed 
by MacGibbon and Weber in ref. fragment 
into hadrons, decaying in their turn into stable 
particles. We have modeled the interactions of 
the quarks and gluons before hadronization takes 
place. In this regime, as suggested by Heck- 
ler in ref. |8), a quark-gluon plasma similar to 
the electron-photon photosphere in QED may de- 
velop, analogously changing the energy spectrum 
of the particles. 

4.5 Test Particle Method in QCD 

Our treatment does not attempt to give a detailed 
model of QCD interactions after hadronization be- 
gins. However, the onset of the photosphere can 
be established in terms of free quarks and glu- 
ons interacting with each other. We assume that 



hadronization occurs at a distance of 



A 



QCD 



and that no significant softening of the particle 
spectrum occurs after this point. Hence, we make 
the same test particle simulation as in QED, only 
with different interaction cross section, masses, 
and number densities. Due to the much larger 
coupling constant and greater number of degrees 
of freedom, we expect the photosphere to develop 
at temperatures Tbh f ar below the critical tem- 
perature for QED, T c ~ 50 GeV, and to reach 
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higher densities than in the QED case. 

To investigate the photosphere in QCD, we 
recall that the key ingredient was the inclusion 
of 2 — > 3 body interactions in the collision term 
of the Boltzmann equation (^). Like electrons 
and photons, quarks and gluons can also inter- 
act via bremsstrahlung (qq — ► qq g) and pair pro- 
duction {qg — > qqq). The dominant diagrams 
(which give large logarithms in the cross section 
at low momentum transfer due to t-channel prop- 
agators) are the same as for QED, fig. 0. In ad- 
dition, there are intrinsically nonabelian contri- 
butions like fig. |l^ which we do not expect to 




Figure 13: Nonabelian contribution to quark- 
gluon bremsstrahlung . 

be parametrically larger than those of fig. |[ We 
will therefore make the simple approximation of 
modeling the bremsstrahlung and pair production 
cross sections using the same form as in QED, but 
replacing the fine structure constant by a s and the 
electron mass by the mass of the lightest quark: 



a 



QCD 
brem 



8a 3 s , 2E 
— ^In . 

m th m th 



(37) 



m Q + m p using the QCD plasma 



and we take m q to be the aver- 



Here mth - y 
mass, eq. (j2j 
age of the up and down quark masses, m q = 8 
McV. Based on the previous QED results, heav- 
ier quarks are expected to make a subdomi- 
nant contribution to the photosphere, since their 
bremsstrahlung cross sections are smaller by a fac- 
tor of the mass ratio squared. The dependence of 
the thermal quark mass on radial distance in the 
photosphere was already shown in fig. |[ 

In contrast to QED, the QCD coupling con- 
stant depends strongly on energy. To leading or- 
der in perturbation theory |l7j ], 
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(33 - 271/) In (jp/A 



2 

QCD 



(38) 



where in our radial evolution we take [i to be the 
average particle energy at a given radius, Aqcd ~ 



200 MeV, and nj ranges between 3 and 5 for the 
BH temperatures we are considering, depending 
on the number of quark species with masses less 
than the average energy at a given distance from 
the horizon. 

When the effective temperature of the quark- 
gluon plasma becomes of order Kqcd, perturba- 
tion theory in a s breaks down, and the quarks 
and gluons hadronize. This process will be dis- 
cussed in Section 4.7. At distances greater than 
Aq CD , instead of a gas of quarks and gluons, we 
will have a plasma of pions and nucleons, which 
can in principle continue to cool through pion 
bremsstrahlung (nn — > nmr or tttt — > ttttt:). How- 
ever the residual strong interactions of pions and 
nucleons are screened, relative to those of the 
quarks and gluons, by confinement. Also the rel- 
evant scale for the range of pion bremsstrahlung 
is m^ 1 , which is much shorter than the range of 
quark-gluon brehmsstrahlung. We expect these 
two reductions in the strength of 2 — » 3 scattering 
processes to make the hadrons much less effec- 
tive than quarks and gluons in perpetuating the 
photosphere. Thus one might anticipate that the 
QCD photosphere ends rather abruptly beyond 
distances of order Aq Cd . 

To semiquantitatively investigate this post- 
hadronization regime, we modeled the hadron gas 
by replacing a s by the pion-nucleon fine-structure 
constant, / 2 /47r = 0.09, obtained from low-energy 
pion-nucleon scattering [Q; we also substituted 
the quark mass with the pion thermal mass, which 
we estimate analogously to eqs. ( |25| ) and (p7[). 
The result is that hadron-hadron interactions are 
indeed ineffectual for prolonging the photosphere. 
Henceforth we simply use the hadronization crite- 
rion to determine where the photosphere ends. 



4.6 Parameters of QCD 
sphere 



Photo- 



We have found that the QCD photosphere starts 
to develop for any black hole whose temperature 
exceeds a critical value 



T r ~ 175 MeV. 



(39) 



This is more than two orders of magnitude lower 
than the critical temperature for the QED photo- 
sphere. It agrees with the analytical estimate in 
ffl, T c ~ Kqcd- Recall that T c is defined to be 
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the temperature at which each particle on average 
undergoes one scattering during its outward prop- 
agation. In Fig. [l4| we show the average number of 
scatterings per particle as a function of distance, 
Af(r), at several BH temperatures. In contrast to 
the QED case (fig. |B|), where Af(r) levels off at 
a universal value of the final radius rj marking 
the end of the photosphere, in QCD rj depends 
strongly on the temperature. The radius at which 
hadronization takes place grows with Tbh , which 
will be quantified below. 
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Figure 14: Average number of scatterings N in 
QCD photosphere as a function of radius r for 
several black hole temperatures. 

The precise definition of the final radius, r/, 
differs in the QCD case from that which we used 
for QED: for QCD we must decide just where 
hadronization occurs. One can imagine several 
possible criteria: when the interparticle spacing 
faO"/)] 3 begins to exceed Aq^; when the 
average particle energy E(rf) becomes of order 
Aqcd; or when the coupling constant a s (rf) be- 
comes of order unity. We find that all three are 
roughly equivalent in that the value of rj depends 
only marginally on which one is used. We adopt 
the coupling constant criterion to define r j in the 
results that follow. Fig. [l5| shows that rj is well 
described by a logarithmic growth in the photo- 
sphere radius with the black hole temperature: 
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1 GeV/ 



(40) 



1.45 ± 0.06 GeV" 



0.29 fm. 



10.0 p 

9.0 - 

8.0 - 

7.0 - 

6.0 - 

5.0 - 

4.0 - 

3.0 - 

2.0 - 

1.0 - 

0.0 - 



— Logarithmic regression 




1 10 
In (T/GeV) 



100 



Figure 15: Radius of the outer surface of the QCD 
photosphere (rf) versus logarithm of the black 
hole temperature. 

QCD also differs from QED for the onset of 
the photosphere, whose inner radius is denoted by 
fj. We find that significant QCD interactions be- 
gin quite close to the horizon: fj ~ for the 
whole range of BH temperatures for which the 
QCD photosphere forms. 

The parameter which best characterizes the 
intensity of interactions in the QCD photosphere 
is the total particle production factor, P(Tbh) = 
N(r f )/N(r h ). As shown in Fig. |l|, we find that 
it increases quite linearly with black hole temper- 
ature: 



P(Z 



BH) 



.62 ±0.01) 



T B 



H 



GeV 



(41) 



This, again, contrasts with the QED case, which 
displayed a noticeable quadratic dependence on 
Tbh at high temperatures. 

The final spectra at the photosphere surface 
are not exactly thermal, but they can be fitted 
over most of the range where dN/dE is significant 
by a Boltzmann distribution, 



dN 
~dE 



oc exp(-£/T ), 



(42) 



where the effective temperature at the photo- 
sphere is independent of Tbh'- 
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Figure 16: Total particle production in the QCD 
photosphere versus black hole temperature. 



T = 300 MeV. 



(43) 



In fig. 17 we show the superimposed spectra at 
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black hole that has reached a temperature greater 
than T c = 175 MeV, corresponding to a mass 
M < 5 x 10 14 g, emits quarks and gluons which 
almost immediately begin interacting to form a 
photosphere very close to the horizon, ^$0.1 
fm. As the black hole temperature continues to 
rise, the photosphere inner radius shrinks along 
with the horizon (r^ ~ = 1/AttTbh), while 
the outer radius grows logarithmically with Tbh- 
Particles emitted from the horizon with average 
energy Ei ~ 3Tbh fragment in the photosphere 
into lower energy particles. The higher Tbh, the 
more particles are created (Fig. |l^). The aver- 
age particle energy decreases as they propagate 
outward, until it reaches Ef ~ 300 MeV at the 
photosphere edge, where hadronization occurs. 

The results for the photosphere parameters 
for several characteristic black hole temperatures 
are summarized in Table ||. As is evident, the av- 
erage particle energy can decrease by several or- 
ders of magnitude in a QCD photosphere. How 
the full spectrum changes is illustrated in fig. 18, 
where we show the energy distributions of the par- 
ticles at the horizon and at the photosphere edge 
of a 1.5 GeV black hole. In this example the av- 
erage energy decreases by a factor of 13, and the 
number of particles increases by the same factor, 
given by P(1.5 GeV) [eq. @)]. 

Of course it is not the quarks and gluons that 
might reach a distant observer, but rather the sub- 
sequently produced hadrons and their decay prod- 
ucts, especially the photons. We will therefore 
consider the processes by which the QCD plasma 
creates a potentially observable gamma ray or an- 
tiproton signal. 



Figure 17: Linear-log plot of QCD photosphere 
spectra (normalized to have same peak value) for 
three different horizon temperatures, and the an- 
alytic fit dN/dE = exp(-£7300 MeV). 

the photosphere for three different initial temper- 
atures, Tbh = 1-5, 15 and 30 GeV, along with the 
fit (Q). The spectra are normalized to have the 
same maximum value so the shapes can be com- 
pared. They rise very sharply from dN/dE = 
at E = 0. The actual normalization of the flux 
dN/dEdt will be discussed below. 

From our simulations at different tempera- 
tures we can reconstruct the time evolution of 
the QCD photosphere of an individual BH. A 



5 Possible Experimental 
Consequences of the 
Photosphere 

Finally we would like to see what difference the 
photosphere makes for observational cosmology or 
astrophysics. From the foregoing it is clear that 
the spectrum of radiation from an individual black 
hole near the end of its existence will be consid- 
erably softened from the usual expectation based 
on the Hawking flux, eq. ([!]). In addition, the 
integrated contributions of black holes to the dif- 
fuse gamma ray or cosmic antiproton backgrounds 
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Ta r*j Tu 


T * 
' J 


J - 1 % 


^ i 


200 MeV 


2.2 fm 


4.9 fm 


600 MeV 


300 MeV 


1 GeV 


0.4 fm 


19 fm 


3.0 GeV 


300 MeV 


50 GeV 


0.008 fm 


49 fm 


156 GeV 


300 MeV 



Table 2: QCD photosphere parameters for several black hole temperatures obtained from the test 
particle simulation. 
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Figure 18: Energy distribution of quarks and glu- 
ons at the horizon of a 1.5 GeV black hole and at 
the edge of the QCD photosphere 



5.1 Hadronization and Final Spec- 
trum of QCD photosphere 

Roughly speaking, the QCD interaction is pertur- 
bative (a s < 1) when the distance between the 
particles is smaller than Ag^ D . This condition is 
satisfied in the photosphere region. At larger dis- 
tances, however, vacuum fragmentation of quarks 
and gluons will occur, which is what happens at 
the photosphere edge. For an accurate calculation 
of the final photon spectrum, we should first com- 
pute the neutral pion flux coming from the photo- 
sphere parton distributions using a jet fragmenta- 
tion code, and then the photon flux resulting from 
7T° — » 77 decays. However, following ref. |ll| we 
can estimate this spectrum as a convolution of 
the quark-gluon spectrum, available from our test 
particle simulation, with the pion fragmentation 
function |^] and the Lorentz-transformed spec- 
trum of photons from 7r° decay :^] 



d 2 N 7 
dE 1 dt 



dEn 



dg„ 7 (E„) d 2 N n 



En 



dE^ 



dErr dt ' 



(44) 



where Eq = E 1 + mt K jAE 1 . The number of pho- 
tons of energy E 7 created by a pion moving with 
velocity /3 and decaying isotropically in its rest 
frame is 



should be shifted to lower energies. The first step 
is to compute how quarks and gluons in the QCD 
photosphere fragment into observable particles. 
We will then integrate the individual BH fluxes 
over time and over the initial mass distribution of 
BH's to arrive at the diffuse background fluxes. 



dg^E v ) = 2 = 2 
dE~f jmirfi y/El - ml ' 



3 We differ with Jis|] concerning the limits of integration 
here. The minimum and maximum photon energies from 
a decaying pion boosted to energy = ^m-^ are E 1 = 

£V(1± yl - 7~ 2 )/2. Inverting this equation implies that 
a photon of energy Ej can come from pions with energies 
satisfying En > E 1 + /4E 7 . 
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where 7 = (1 — /3 2 ) 1//2 . The pion flux is 
d 2 N, 



dE-x dt 



J , n dg ll (Q 1 E^)fN i 
^ dE^ dQdt ' [ Dj 



where the sum is over relevant species in the 
plasma (quarks and gluons). For the number of 
pions with energy E^ produced per unit energy 
by each quark or gluon, we use the empirical frag- 
mentation function lEU 01 



dg 3V (Q,E v ) 
dEn 




Et, 



(47) 
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photosphere 
Excluding 
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Lastly, d 2 Nj/dQdt is the quark or gluon flux at 
the outer edge of the photosphere. In what follows 
we treat quark and gluon jets equally and thus 
write V\ d 2 Nj/dQdt = dN/dQ. 

Combining (|44|) through (f47j ) we obtain: 

d 2 N~ f°° , 15/8 

dE 1 dt J Eo El /2 ^jEl-ml 

where E Q = E 1 + m\ / AE 1 . 

We have calculated the integral (ji|) numer- 
ically for a large range of black hole tempera- 
tures. The results for one of them (Tbh — 50 
GeV) are presented in Fig. |l^ and compared to 
the results obtained neglecting the photosphere, 
but taking into account direct quark fragmenta- 
tion at the horizon and subsequent tt° decay, as 
in ref. g. Also shown are the spectra of photons 
emitted directly from the black hole (neglecting 
the QED photosphere) and from the QED photo- 
sphere, which just starts to form at this tempera- 
ture. The actual full spectrum of a 50 GeV black 
hole is the addition of the two solid lines. The re- 
sults are in agreement with fig. 1 of ref. jl9) , except 
for the QED photosphere spectrum, as discussed 
above. 

We see that the photon flux from the QCD 
photosphere peaks at an energy of if the 

pions were at rest then all photons would have 
this energy, but E y is Doppler shifted by the pion 
motion, and the width of the distribution grows 
with the average energy of the pion. Figure ^0] 
shows the photon spectra for a large range of BH 



Figure 19: Photon emission spectrum from T = 
50 GeV (M = 2 x 10 12 g) black hole. Solid lines are 
spectra which include photospheres. Dashed lines 
are given for comparison and represent the results 
for photons from direct quark fragmentation at 
the horizon and subsequent ir° decay (QCD), and 
for direct photon emission neglecting the QED 
photosphere (QED). 




Figure 20: Photon emission spectrum from ir° de- 
cay for black holes, ignoring photosphere, with 
Tbh = 50 MeV (narrowest), increasing Tbh by 
factors of 10, up to 500 TeV (widest). Dashed 
curve shows the spectrum including photosphere 
(see text for normalization). 

temperatures in the case where the photosphere 
is neglected. As a function of y = log(2E 7 /m 7T ), 
these spectra are symmetric under y —> —y. Be- 
cause the total power output of the BH goes like 
/ dQ QdN/dQ ~ T BH , it is convenient to display 
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instead the quantity Tg H dN~,/dEj, as we have 
done in Fig. |2(]. We fit these functions to 10th 
order polynomials, 

log 10 (T^dN y /dE^ ^J2 a n\v\ n l 

i 

y = log 10 2E 1 /m 7 t (49) 

whose coefficients (up to 8th order, which is still 
a good approximation) are given in table |. The 
absolute normalization follows from the ansatz 
dN/dQ = (Q/T BH ) 2 exp(-Q/T Bff )/(3!) which 
we took as an approximation to the actual Hawk- 
ing distribution (up to an overall normalization 
factor which can be computed) to obtain these 
functions. They will be useful later when we in- 
tegrate the BH radiation over time. Although the 
spectra are really analytic at y = 0, the deriva- 
tive changes so quickly there that the best fit is 
obtained by using an odd function for y > and 
then letting y — > \y\ to cover the y < side. Also 
shown in fig. |2^ is the corresponding spectrum, 
divided by Tq , for the photosphere, whose parton 
flux is taken to be dN/dQ = exp(— Q/Tq), where 
To = 300 MeV. Like the nonphotosphere results, 
a photon spectrum computed using this distribu- 
tion should be multiplied by a factor of T BH (for 
Tbh > Tq) to represent a BH with horizon tem- 
perature Tbh, to insure that the power output at 
the photosphere is equal to that at the horizon for 
any value of Tbh- 



Finding the contribution of black hole radiation 
to the diffuse photon background consists of two 
steps: (1) first integrate the contribution of a sin- 
gle BH over time, and (2) integrate this result 
over the initial mass distribution of the individual 
BH's. Let us denote the time-integrated spectrum 
emitted by a single black hole by dN^/dE 1 . Go- 
ing from an initial time ti to the final time i/, and 
accounting for the redshifting of the photons be- 
tween the time of emission and the present (to), 



integrate over the BH temperature, however. The 
time-temperature relation can be found by equat- 
ing the rate of change of the BH mass with its 
power output, to obtain EL || 



dT 

~dt 

a 



a(T)GT 4 , 

0.57 di /2 + 0.23 tfi, 



(51) 



where d s is the number of degrees of freedom 
(spin, charge and color) of spin s which can be 
emitted by the BH at the given temperature. If 
we ignore the weak T-dependence of a, eq. (|51 
gives the time-temperature relation 



T- & ) 



(52) 



for a BH with initial temperature Ti at time ti. 
Let us now assume that ti = and define T* as 
the initial temperature of a black hole which is 
disappearing today, T* « 100 MeV. Then <o = 
l/(3aT%), and we can write the redshift factor as 



Z = 1+z 



2/3 




Thus 

dNi f Tf{Ti) dT 



5.2 Diffuse gamma ray background dE^ 



Z(T)^[Z(T)E y }. (54) 



aGT 4 v ' dE. 



The final temperature in this expression is given 

by 



Tf 



Ti(l - (T,/T,) 3 )- 1 /3 j Ti < n 
oo, Ti > T* 



(55) 



Next we must add up the contributions from 
all black holes. The distribution of initial BH 
masses is taken to be [0] 



dNi 
dTL 



"f dN 
dtZ(t)—t[Z(t)E], 



dE^ 



(50) dN _ {13 - 2)n BH Pc 



where dN 7 /dE 7 is the 



77 flux derived in 



the previous section. Notice that we must multi- 
ply both dN~[/dEj and its argument by the red- 



dMi 



Ml 



Mi 
Ml 



= C M M U 



(56) 



where VL B h is the fraction of the critical density 
of the universe which is in primordial black holes, 



shift factor Z = (1 + z). It is more convenient to (3 — 2.5 for the usual equation of state p = p/3, 
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Tbh 


ao 




0,2 


03 


04 


05 


a,Q 


07 


a$ 


50 MeV 


-5.660 


-1.634 


-3.055 


-1.853 


1.308 


-1.518 


-1.207 


2.169 


-1.596 


500 MeV 


-4.938 


-0.973 


-1.059 


-0.353 


1.544 


-2.191 


1.566 


-0.710 


0.189 


5 GeV 


-5.289 


-0.822 


-1.124 


1.016 


-0.818 


0.6191 


-0.418 


0.195 


-0.057 


50 GcV 


-5.882 


0.237 


-4.838 


7.965 


-8.251 


5.528 


-2.430 


0.695 


-0.125 


500 GcV 


-6.270 


-0.795 


-1.136 


1.105 


-0.747 


0.378 


-0.146 


0.041 


-0.008 


5 TcV 


-6.770 


-0.777 


-1.243 


1.338 


-0.991 


0.514 


-0.184 


0.044 


-0.007 


50 TcV 


-7.270 


-0.735 


-1.402 


1.549 


-1.105 


0.521 


-0.163 


0.033 


-0.004 


500 TeV 


-7.770 


7.316 


-24.23 


27.72 


-17.30 


6.561 


-1.580 


0.244 


-0.023 


photosphere 


-5.070 


-1.177 


-1.629 


-0.515 


1.506 


-1.965 


1.155 


-0.411 


0.0648 



Table 3: Coefficients for photon fluxes from ir° — ► 77 used in eq. (|49|). 



and M* » 10 15 g is the mass of a BH with T. t = T* 
(hence a BH which is disappearing in the present 
epoch). Then, trading Mj for Tj, the integral over 
initial BH's gives the spectrum of diffuse gamma 
rays as 



dN-y 

~d~EL 



= C 



A I 



dT, T- 



j3-2 



dN x 

IeZ' 



(57) 



The photon flux per unit energy is ci/V 7 /d£ , 7 
times the speed of light and a geometric factor, 
(471-)- 1 d(j> J,, 1 d(cos 6>) cos 6 = 1 /4. We have 
computed this flux both with and without the 
QCD photosphere to see the effect of the latter. 
The result is shown in fig. [2l], where we have nor- 
malized the no-photosphere curve to agree with 
the predictions of ref. j§ at E = 100 MeV, for the 
case of Qbh = 7.6 x 10~ 9 /iq 2 , which saturates 
the experimental limit. One sees that although 
the photosphere dramatically suppresses the spec- 
trum for E > 100 MeV, the effect is small at lower 
energies. 

The 100 MeV energy range is the most im- 
portant one for setting limits on the primordial 
BH contribution to the energy density of the uni- 
verse, for the following reason. The theoretical 
prediction for the diffuse background spectrum 
goes like E~ x at low energies and £V 3 at high 
energies, with E ~ 100 MeV ~ m, being the re- 
gion where the slope changes. On the other hand, 
the extragalactic flux measured by the Energetic 
Gamma Ray Experiment Telescope (EGRET) has 
an E^ 2 - 1±om energy spectrum |^|, intermediate 
between the two theoretical slopes. Therefore as 
one increases £Ibh from small values, the first 
place where the theoretical prediction comes into 
conflict with the observation is at the "knee" of 



> 

CD 



CO 



o 
o 



o> ■ 
o 




0.5 1.0 1.5 2.0 2.5 3.0 3.5 4.i 
Log (Photon Energy / MeV) 



Figure 21: Theoretical prediction for the diffuse 
gamma ray spectrum from decaying black holes, 
with and without the QCD photosphere, assuming 
n BH = 7.6 x lQ- 9 hn 2 . 



the theoretical spectrum. The photosphere has 
only a 60% effect on the flux at these energies. 
In figure ^ we show where a line with the ob- 
served slope first encounters the predicted spec- 
tra as Qbh is increased. The intercept decreases 
by a factor of 10 2 = 1.6 when the photosphere 
is taken into account; thus the bound on CIbh is 
only slightly weakened. 

5.3 Antiproton background 

A similar effect of the photosphere can be found 
in the predicted flux of antiprotons, which is inter- 
esting because of two current experiments which 
search for antimatter coming from cosmic sources, 
BESS II and AMS ||. Since the effective tem- 
perature of the BH is limited to ~ 300 MeV by 
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Eq. (58) only gives a good approximation for val- 
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Figure 22: Close-up of the 100 MeV region in 
the previous figure, showing where a line with the 
slope of the observed spectrum becomes tangent 
to each curve. 



the QCD photosphere, we expect a reduction in 
the flux of protons and antiprotons relative to pre- 
dictions using the Hawking spectrum. 

The computation is somewhat simpler than 
for photons since we only need the fragmentation 
of quark and gluon jets into antiprotons, with 
no additional subequent decay as in the case of 
7T° — * 77. A rough fit to the fragmentation func- 
tion can be inferred from actual data for e + e~ — > 
jets — » hadrons. Let Q be the quark or gluon en- 
ergy, p p the antiproton momentum, and x p = p/Q 
the momentum fraction. From fig. 17 of ref. |p5[ 
we find that the normalized cross section for p, p 
production can be fit by 



da 



mi In x p + bi , x p > x T 



a tot dx p \ m 2 In x p + b 2 , x p < x 

T. 



where 



mi = -0.259, 61 = 0.014; 

m 2 = 0.318, 62 — (m-2 -niijlni, 



(58) 



(59) 
+ 61, 



and the momentum fraction where the distribu- 
tion peaks, x m ax, empirically depends on the par- 
ton energy according to 



In x„ 



-0.97 log(Q/GeV) - 0.95. (60) 



ues such that (x p /atot)da/dx p £ 0.03; the tails of 
the distribution are better represented by a gaus- 
sian, 



x p da 
a tot dx p 



F(Q) exp 



In 2 x r 



c / 



° 2 {Q) 



(61) 



where the width is supposed to depend on Q 
like a 2 {Q) = C(ln 3/2 4Q 2 /A 2 - ln 3/2 fi 2 /A 2 ), with 
fi S 0.35 GeV, A = A QCD = 0.2 GeV, F(Q) 
chosen so that J Q dx p (x p /atot)da/dx p = 1 and C 
being a constant. However we did not find this 
to be a good representation of the actual data in 
the vicinity of the peak, for any constant value 
of C. We have thus relegated the form (|6l|) for 
representing the tails of the distribution, which in 
any case make a subdominant contribution to the 
final antiproton flux. 

The fragmentation function is related to the 
cross section by 



dg 



1 da 



E Xp da 



dEp a tot dE p p 2 a tot dx p 
The p, p flux from a single black hole is then 



(62) 



dNp 
dE n 



do , „ „ „ dN 



dE n 



dQ 



(63) 



We show the numerical results for four different 
cases in fig. [2^: black holes with T = 1 GeV and 
T = 10 GeV, ignoring the photosphere, and the 
same BH's taking the photosphere into account. 
We have used the results shown in fig. 1 of ref. 
Q to fix the absolute normalization. The only 
effect of the horizon temperature, when the pho- 
tosphere effects are included, is to multiply the 
distribution by a factor of (T/T ) 2 for T > T , 
where To = 300 MeV is the effective temperature 
of the photosphere. This factor comes from de- 
manding that the total power output of the BH is 
the same with or without the photosphere. 

To find the diffuse antiproton flux, we should 
integrate over the black hole temperature and the 
distribution of initial masses (temperatures) as we 
did for photons, eqs. (|54|) and (|57|). The only dif- 
ference is that we are interested in nonrelativistic 
as well as relativistic protons, so we must redshift 
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Figure 23: Instantaneous p,p fluxes from BH's 
with T B h = 1 GeV (solid lines) and 10 GeV 
(dashed lines), including ("w/ph.") or neglect- 
ing ("no. ph.") the QCD photosphere, plotted as 
a function of the kinetic energy. 



the momentum rather than the energy. Instead of 
the factor ZdN 1 / dE^[ZEi\ in eq. (fl), we get 



dN v . 

z ii} ZE ^ 



Z 2 E n dN, 



E dE p 



(64) 



where Eq = y Z 2 (E 2 — m 2 ) + m 2 . The result 

(normalized to agree with fig. 2 of ref. Q) is shown 
in fig. ^4|. Again, we take the BH density to be the 
maximum allowed by the gamma ray background, 
n BH = 7.6 x 1CT 9 V 2 - 

The recent observations by BESS give an an- 
tiproton flux of (8±2) x 10~ 3 m~ 2 sr^scc^GeV -1 
at E p = 0.2 — 0.3 GeV, and higher values up to 
2 x 10~ 2 at larger energies; thus the bound on 
flsH from antiprotons is set by the low energy 
range of the observations. Comparison with fig. 
^l| shows that the nonphotosphere prediction for 
the p flux is somewhat in excess of the data at 
these energies, suggesting that one might be able 
to set slightly stronger limits on Qbh using the an- 
tiproton flux rather than gamma rays. However, 
the predicted flux taking into account the photo- 
sphere is 10°' 6 = 4 times weaker, and in better 
agreement with the data. 

It might seem surprising that the low energy 
p flux is not degraded more than it is by the pho- 



Figure 24: Diffuse background p,p flux from inte- 
grating over time and initial distribution of BH's, 
as a function of the kinetic energy. Qbh is the 
same as in fig. CU 



tosphere at kinetic energies below m p . Apparently 
it is the tails of the distributions of the underlying 
partons which are mostly responsible for produc- 
ing low-energy protons, and the largest contribu- 
tion will come from black holes with temperatures 
somewhat below the nucleon mass. In this regime 
the difference between having the QCD photo- 
sphere or not is minimized. 



6 Conclusion 

Our main results can be summarized as follows. 
The test-particle method of solving the Boltz- 
mann equation, previously used for analyzing 
heavy ion collisions, was applied to the problem of 
black hole evaporation. The method was adapted 
to the situation of steady-state diffusion of parti- 
cles emitted by a microscopic black hole. A code 
to simulate the bremsstrahlung and pair produc- 
tion interactions of the test particles was devel- 
oped, leading to solutions for the particle distri- 
bution functions at any distance from the black 
hole horizon. 

Simulation of microscopic black hole emission 
in both QED and QCD energy ranges corrobo- 
rates the idea of photosphere formation pioneered 
by Heckler in We find that any black hole of 
mass M < 5 X 10 14 g develops a cloud of interact- 
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ing quarks and gluons which extends a certain dis- 
tance from the black hole horizon. The evolution 
of such small (r < 0.08 fm) black holes is dom- 
inated by mass loss through Hawking radiation. 
Part of this radiation is in free quarks and gluons 
which are processed in the QCD photosphere un- 
til their average energy drops to the point E ~ 
h-QCD, where they hadronize into stable particles 
and fast-decaying pions. Another part consists of 
electrons, positrons and photons. Once the black 
hole mass drops below M ~ 2 x 10 12 g, these 
particles interact significantly enough to form an- 
other, less dense cloud at a distance about 700 
times the horizon radius. This QED photosphere 
extends over a distance of about 400 fm, where it 
dissipates and emits much less energetic, but more 
numerous, electrons, positrons and photons. 

Energy distributions of the particles leav- 
ing both photospheres were obtained and shown 
to greatly differ from the original nearly-thermal 
Hawking distributions by being softened to much 
lower average energies: E ~ 300 MeV for QCD, 
E ranging from 100 GeV for a 10 12 g black hole 
- 0.5 MeV for a 10 6 g black hole in the QED case. 
We used the QCD spectra to compute the con- 
tributions of individual black holes and all BH's 
in the universe to potentially observable gamma 
ray and antiproton signals, and compared to the 
previous expectations based on ignoring the pho- 
tosphere. In the regions where the experimental 
sensitivity is greatest, the photosphere lowers the 
fluxes by only a small factor: 1.6 for photons and 
4 for antiprotons. 

Our findings do not support the approxima- 
tion made in ref. || of treating the photosphere as 
a fluid. Rather we get a picture of a steadily ex- 
panding cloud of particles which never quite ther- 
malizes, and has interactions which are primarily 
low in momentum transfer. This is how we in- 
terpret the fact that our Boltzmann code gives 
much smaller photosphere, hence much less en- 
ergy degradation of particles, in the case of the 
QED photosphere, than claimed in J^]. This dis- 
crepancy did not appear for the QCD case be- 
cause there both approaches find that the photo- 
sphere ends when hadronization begins. Further- 
more we do not find a relativistically expanding 
photosphere, which was also claimed in || on the 
basis of the fluid approach. This led us to find 
larger differences in the diffuse gamma ray back- 
ground between the photosphere and nonphoto- 



sphere predictions than found by ref. fi9f . The 
reason is that |l9| undoes the energy-degrading 
effects of the photosphere to a large extent by 
boosting the distributions from the fluid frame to 
the observer frame, a step which is not necessary 
in our approach, since we always work in the latter 
reference frame. 

It is disappointing that the observable con- 
sequences of the photosphere are small in the ex- 
perimentally interesting energy ranges. If we were 
lucky enough to have a nearby BH reach the end 
of its existence however, a real experimental test 
might be possible, since the spectra for individual 
BH's have radically different characteristics with 
or without a photosphere. On the theoretical side, 
the question of whether black hole photospheres 
indeed develop hinges crucially on whether the 
range of the bremsstrahlung interaction is really 
l/m e or something effectively shorter. Although 
the arguments supporting this claim look plausi- 
ble, it is perhaps deserving of more detailed study. 
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